Last updated: 2024-07-17

Checks: 6 1

Knit directory: CocaineIVSA_RNASequencing_mDANeurons/

This reproducible R Markdown analysis was created with workflowr (version 1.7.0). The Checks tab describes the reproducibility checks that were applied when the results were created. The Past versions tab lists the development history.


Great! Since the R Markdown file has been committed to the Git repository, you know the exact version of the code that produced these results.

Great job! The global environment was empty. Objects defined in the global environment can affect the analysis in your R Markdown file in unknown ways. For reproduciblity it’s best to always run the code in an empty environment.

The command set.seed(20240716) was run prior to running the code in the R Markdown file. Setting a seed ensures that any results that rely on randomness, e.g. subsampling or permutations, are reproducible.

Great job! Recording the operating system, R version, and package versions is critical for reproducibility.

Nice! There were no cached chunks for this analysis, so you can be confident that you successfully produced the results during this run.

Using absolute paths to the files within your workflowr project makes it difficult for you and others to run your code on a different machine. Change the absolute path(s) below to the suggested relative path(s) to make your code more reproducible.

absolute relative
/home/avm27/Documents/Raw_Sequencing_Data/CIVSA2/CocaineIVSA_RNASequencing_mDANeurons .

Great! You are using Git for version control. Tracking code development and connecting the code version to the results is critical for reproducibility.

The results in this page were generated with repository version ca44efa. See the Past versions tab to see a history of the changes made to the R Markdown and HTML files.

Note that you need to be careful to ensure that all relevant files for the analysis have been committed to Git prior to generating the results (you can use wflow_publish or wflow_git_commit). workflowr only checks the R Markdown file, but you know if there are other scripts or data files that it depends on. Below is the status of the Git repository when the results were generated:


Ignored files:
    Ignored:    .Rproj.user/

Untracked files:
    Untracked:  VennDiagram.2024-07-17_07-40-38.log
    Untracked:  VennDiagram.2024-07-17_07-40-39.log
    Untracked:  VennDiagram.2024-07-17_07-58-38.log
    Untracked:  VennDiagram.2024-07-17_07-58-44.log
    Untracked:  VennDiagram.2024-07-17_07-59-00.log
    Untracked:  VennDiagram.2024-07-17_08-04-55.log
    Untracked:  VennDiagram.2024-07-17_08-11-08.log
    Untracked:  data/CIVSA_normalized_counts_DESeq2.txt
    Untracked:  data/gene_count_matrix_CIVSA.csv
    Untracked:  data/sample_info_CIVSA.csv
    Untracked:  output/Sig_Res_Craving_vs_FoodTrained.csv
    Untracked:  output/Sig_Res_Craving_vs_Maintenance.csv
    Untracked:  output/Sig_Res_Maintenance_vs_FoodTrained.csv
    Untracked:  plots/
    Untracked:  results/

Unstaged changes:
    Modified:   analysis/_site.yml

Note that any generated files, e.g. HTML, png, CSS, etc., are not included in this status report because it is ok for generated content to have uncommitted changes.


These are the previous versions of the repository in which changes were made to the R Markdown (analysis/index.Rmd) and HTML (docs/index.html) files. If you’ve configured a remote Git repository (see ?wflow_git_remote), click on the hyperlinks in the table below to view the files as they were in that past version.

File Version Author Date Message
Rmd ca44efa avm27 2024-07-17 Testing
html 21b70e3 avm27 2024-07-17 Build site.
Rmd ddb0db7 avm27 2024-07-17 Testing
html a1c090c avm27 2024-07-17 Build site.
Rmd adbd778 avm27 2024-07-17 Testing
html dc04b15 avm27 2024-07-17 Build site.
Rmd 65a97a6 avm27 2024-07-17 Testing
html e58ce5e avm27 2024-07-16 Build site.
Rmd ad8ae96 avm27 2024-07-16 Start workflowr project.

Microglia contribute to methamphetamine reinforcement and reflect persistent transcriptional and morphological adaptations to the drug

Samara J. Vilca1,2†, Alexander V. Margetts1-3†, Isabella Fleites1-3, Claes R. Wahlestedt1-3 & Luis M. Tuesta1-3*

1 Department of Psychiatry & Behavioral Sciences 2 Center for Therapeutic Innovation 3 Sylvester Comprehensive Cancer Center University of Miami Miller School of Medicine, Miami, FL 33136 † Denotes equal contribution * Corresponding author ()

1 Abstract

Abstract Methamphetamine use disorder (MUD) is a chronic, relapsing disease that is characterized by repeated drug use despite negative consequences for which there are currently no FDA approved cessation therapeutics. Repeated methamphetamine (METH) use induces long-term gene expression changes in brain regions associated with reward processing and drug-seeking behavior, and recent evidence suggests that methamphetamine-induced neuroinflammation may also shape behavioral and molecular responses to the drug. Microglia, the resident immune cells in the brain, are principal drivers of neuroinflammatory responses and contribute to the pathophysiology of substance use disorders. Here, we investigated transcriptional and morphological changes in striatal microglia in response to methamphetamine-taking and during methamphetamine abstinence, as well as their functional contribution to drug-taking behavior. We show that methamphetamine self-administration induces transcriptional changes related to protein folding, mRNA processing, immune signaling, and neurotransmission in striatal microglia. Importantly, many of these transcriptional changes persist through abstinence, a finding supported by morphological analysis. Functionally, we report that microglial ablation increases methamphetamine-taking, possibly involving neuroimmune and neurotransmitter regulation. In contrast, microglial depletion did not alter methamphetamine-seeking behavior following 21 days of abstinence, highlighting the complexity of drug-seeking behaviors. Taken together, these results suggest that methamphetamine induces both short and long-term changes in striatal microglia that contribute to altered drug-taking behavior and may be leveraged for preclinical development of methamphetamine cessation therapeutics.

2 DESeq2 Analysis

2.1 Setup

Building DESeq2 input files from gene count matrices output by StringTie following removal of samples that are outliers.

Running DeSeq2 analysis with factor levels of FoodTrained vs Maintenance, Craving vs Maintenance and Craving vs Foodtrained

3 Plot Estimates and Gather Results

3.1 Gather Results

These results include overview for each result file and dispersion estimates based on count values.

3.1.1 Craving vs FoodTraining

Craving vs FoodTrained CIVSA Results Overview
baseMean log2FoldChange lfcSE stat pvalue padj
Mrpl15|Mrpl15 2162.26767 0.3962276 0.3468265 1.1424376 0.2532722 0.9784152
Gm39587|Gm39587 20.30129 -2.1136742 3.3911278 -0.6232954 0.5330904 0.9997970
Lypla1|Lypla1 994.90242 0.0198799 0.8799100 0.0225931 0.9819748 0.9997970
Xkr4|Xkr4 4459.42848 -0.0337912 0.2715766 -0.1244260 0.9009780 0.9997970
Gm39585|Gm39585 5797.38480 0.0861245 0.2414588 0.3566843 0.7213282 0.9997970
Tcea1|Tcea1 2657.79689 -0.4037632 0.3330614 -1.2122786 0.2254058 0.9485680
Rgs20|Rgs20 422.79650 -0.1889366 0.7619761 -0.2479561 0.8041684 0.9997970
Gm16041|Gm16041 21.30106 7.5945314 3.6799019 2.0637864 0.0390380 0.6410945
Atp6v1h|Atp6v1h 10659.10495 -0.4956230 0.2967682 -1.6700675 0.0949060 0.8425404
LOC108167788|LOC108167788 11.33585 -0.8220772 3.4983015 -0.2349933 0.8142140 0.9997970
Npbwr1|Npbwr1 924.87540 -0.3079659 0.6636912 -0.4640199 0.6426335 0.9997970
Oprk1|Oprk1 15661.74234 -0.1442457 0.3231485 -0.4463757 0.6553259 0.9997970
St18|St18 279.69289 -2.0082890 0.9250037 -2.1711146 0.0299225 0.5845690
Rb1cc1|Rb1cc1 12217.86115 -0.0796796 0.2844821 -0.2800866 0.7794110 0.9997970
4732440D04Rik|4732440D04Rik 5369.16278 -0.4918180 0.2427412 -2.0261004 0.0427545 0.6632640
Gm26901|Gm26901 163.46252 -0.2293404 0.9562622 -0.2398300 0.8104620 0.9997970
Pcmtd1|Pcmtd1 4321.05905 0.2034726 0.2645104 0.7692425 0.4417494 0.9997970
Rrs1|Rrs1 625.86245 -0.4850563 0.6605362 -0.7343371 0.4627433 0.9997970
Mybl1|Mybl1 184.22714 1.4451732 1.2066675 1.1976566 0.2310507 0.9561403
Adhfe1|Adhfe1 452.68051 -0.1536903 0.7879557 -0.1950494 0.8453543 0.9997970

out of 23922 with nonzero total read count adjusted p-value < 0.1 LFC > 0 (up) : 250, 1% LFC < 0 (down) : 48, 0.2% outliers [1] : 2272, 9.5% low counts [2] : 1352, 5.7% (mean count < 2) [1] see ‘cooksCutoff’ argument of ?results [2] see ‘independentFiltering’ argument of ?results

3.1.2 Maintenance vs FoodTraining

Craving vs FoodTrained CIVSA Results Overview
baseMean log2FoldChange lfcSE stat pvalue padj
Mrpl15|Mrpl15 2162.26767 0.6589546 0.3707556 1.7773286 0.0755142 0.3782100
Gm39587|Gm39587 20.30129 -1.7619611 3.6258386 -0.4859458 0.6270056 0.9345721
Lypla1|Lypla1 994.90242 0.0914200 0.9407122 0.0971817 0.9225821 0.9949371
Xkr4|Xkr4 4459.42848 0.6341912 0.2902594 2.1849118 0.0288953 0.2180911
Gm39585|Gm39585 5797.38480 1.1798031 0.2580070 4.5727554 0.0000048 0.0001973
Tcea1|Tcea1 2657.79689 0.1319203 0.3560044 0.3705581 0.7109667 0.9583981
Rgs20|Rgs20 422.79650 0.3482376 0.8144256 0.4275868 0.6689520 0.9469053
Gm16041|Gm16041 21.30106 7.1396504 3.9105846 1.8257246 0.0678918 0.3580246
Atp6v1h|Atp6v1h 10659.10495 0.1185235 0.3172387 0.3736100 0.7086945 0.9576220
LOC108167788|LOC108167788 11.33585 -0.5355535 3.7403736 -0.1431818 0.8861466 0.9949371
Npbwr1|Npbwr1 924.87540 -0.6449281 0.7097813 -0.9086292 0.3635459 0.7811879
Oprk1|Oprk1 15661.74234 0.2927281 0.3454541 0.8473720 0.3967878 0.8053124
St18|St18 279.69289 -1.0473877 0.9882176 -1.0598756 0.2892012 0.7120102
Rb1cc1|Rb1cc1 12217.86115 0.3066224 0.3041177 1.0082359 0.3133412 0.7348823
4732440D04Rik|4732440D04Rik 5369.16278 -0.1463996 0.2594995 -0.5641616 0.5726441 0.9132695
Gm26901|Gm26901 163.46252 -0.7535570 1.0236447 -0.7361509 0.4616389 0.8511309
Pcmtd1|Pcmtd1 4321.05905 0.5461501 0.2827538 1.9315391 0.0534164 0.3103856
Rrs1|Rrs1 625.86245 -0.1126443 0.7061219 -0.1595253 0.8732551 0.9949371
Mybl1|Mybl1 184.22714 -2.0915101 1.2989231 -1.6101878 0.1073569 0.4493503
Adhfe1|Adhfe1 452.68051 -1.0496649 0.8432838 -1.2447350 0.2132292 0.6256462

out of 23922 with nonzero total read count adjusted p-value < 0.1 LFC > 0 (up) : 628, 2.6% LFC < 0 (down) : 1046, 4.4% outliers [1] : 2272, 9.5% low counts [2] : 2713, 11% (mean count < 6) [1] see ‘cooksCutoff’ argument of ?results [2] see ‘independentFiltering’ argument of ?results

3.1.3 Craving vs Maintenance

Craving vs Maintenance CIVSA Results Overview
baseMean log2FoldChange lfcSE stat pvalue padj
Mrpl15|Mrpl15 2162.26767 -0.2627261 0.3466449 -0.7579115 0.4485040 0.8651323
Gm39587|Gm39587 20.30129 -0.3516955 3.3976607 -0.1035111 0.9175574 0.9939251
Lypla1|Lypla1 994.90242 -0.0715356 0.8799222 -0.0812976 0.9352053 0.9951721
Xkr4|Xkr4 4459.42848 -0.6679818 0.2714805 -2.4605151 0.0138738 0.1495170
Gm39585|Gm39585 5797.38480 -1.0936780 0.2412810 -4.5327985 0.0000058 0.0002749
Tcea1|Tcea1 2657.79689 -0.5356827 0.3330620 -1.6083574 0.1077569 0.4613839
Rgs20|Rgs20 422.79650 -0.5371710 0.7618243 -0.7051114 0.4807409 0.8833096
Gm16041|Gm16041 21.30106 0.4549475 3.4906011 0.1303350 0.8963014 0.9909408
Atp6v1h|Atp6v1h 10659.10495 -0.6141458 0.2967693 -2.0694383 0.0385050 0.2733236
LOC108167788|LOC108167788 11.33585 -0.2865007 3.5020543 -0.0818093 0.9347983 0.9951721
Npbwr1|Npbwr1 924.87540 0.3369643 0.6640179 0.5074627 0.6118302 0.9372345
Oprk1|Oprk1 15661.74234 -0.4369729 0.3231413 -1.3522657 0.1762903 0.5932565
St18|St18 279.69289 -0.9608973 0.9259208 -1.0377748 0.2993749 0.7560424
Rb1cc1|Rb1cc1 12217.86115 -0.3863013 0.2844705 -1.3579665 0.1744743 0.5910084
4732440D04Rik|4732440D04Rik 5369.16278 -0.3454179 0.2427896 -1.4227046 0.1548218 0.5589467
Gm26901|Gm26901 163.46252 0.5242196 0.9577914 0.5473213 0.5841580 0.9274722
Pcmtd1|Pcmtd1 4321.05905 -0.3426769 0.2644209 -1.2959526 0.1949918 0.6240924
Rrs1|Rrs1 625.86245 -0.3724096 0.6606702 -0.5636845 0.5729689 0.9244236
Mybl1|Mybl1 184.22714 3.5366832 1.2149336 2.9110094 0.0036026 0.0603776
Adhfe1|Adhfe1 452.68051 0.8959769 0.7889411 1.1356702 0.2560946 0.7074992

out of 23922 with nonzero total read count adjusted p-value < 0.1 LFC > 0 (up) : 1084, 4.5% LFC < 0 (down) : 324, 1.4% outliers [1] : 2271, 9.5% low counts [2] : 2713, 11% (mean count < 6) [1] see ‘cooksCutoff’ argument of ?results [2] see ‘independentFiltering’ argument of ?results

3.2 Annotate Results Files and Pull Significant DEGs

This subsection is for annotation the results files to include all necessary information for downstream analyses.

After annotation I subset the results to only view the significantly differentially expressed genes with an adjusted p-value > 0.05.

## Subset for significant genes only
results_sig <- subset(res, padj < 0.05)
results_sig2 <- subset(res2, padj < 0.05)
results_sig3 <- subset(res3, padj < 0.05)

3.2.1 Craving vs FoodTrained

Significant DEGs Craving vs FoodTrained
baseMean log2FoldChange lfcSE stat pvalue padj symbol description entrez ensembl
Gm28836|Gm28836 24.57062 18.512895 3.6825760 5.027159 0.0000005 0.0000748 Gm28836 predicted gene 28836 102635593 NA
Vwa3b|Vwa3b 68.09876 9.322079 2.2726944 4.101774 0.0000410 0.0048096 Vwa3b von Willebrand factor A domain containing 3B 70853 ENSMUSG0….
4930439A04Rik|4930439A04Rik 239.71050 6.086391 1.0543294 5.772760 0.0000000 0.0000025 4930439A04Rik RIKEN cDNA 4930439A04 gene 78119 NA
LOC108167622|LOC108167622 55.54376 21.271643 2.7752955 7.664641 0.0000000 0.0000000 LOC108167622 NA NA
Pard3bos2|Pard3bos2 48.12983 -22.863164 3.6323379 -6.294338 0.0000000 0.0000002 Pard3bos2 par-3 family cell polarity regulator beta, opposite strand 2 105243940 ENSMUSG0….
Col4a3|Col4a3 29.96706 19.601203 3.5629545 5.501390 0.0000000 0.0000079 Col4a3 collagen, type IV, alpha 3 12828 ENSMUSG0….
Dnajb3|Dnajb3 70.36681 9.093523 2.0877638 4.355628 0.0000133 0.0017489 Dnajb3 DnaJ heat shock protein family (Hsp40) member B3 15504 ENSMUSG0….
Asb18|Asb18 12.85015 18.822115 3.6835623 5.109759 0.0000003 0.0000497 Asb18 ankyrin repeat and SOCS box-containing 18 208372 ENSMUSG0….
Gm41920|Gm41920 91.39199 20.766926 3.6195797 5.737386 0.0000000 0.0000030 Gm41920 predicted gene, 41920 105246666 NA
LOC105246710|LOC105246710 41.29964 21.372039 3.6789651 5.809253 0.0000000 0.0000022 LOC105246710 NA NA
LOC108167732|LOC108167732 59.78860 21.982685 3.2366754 6.791749 0.0000000 0.0000000 LOC108167732 NA NA
Fam129a|Fam129a 58.26340 9.542099 2.6179113 3.644928 0.0002675 0.0256134 Fam129a NA NA
Fmo1|Fmo1 51.26689 21.615385 2.7780474 7.780783 0.0000000 0.0000000 Fmo1 flavin containing monooxygenase 1 14261 ENSMUSG0….
Gorab|Gorab 666.84329 2.774594 0.5327144 5.208407 0.0000002 0.0000307 Gorab golgin, RAB6-interacting 98376 ENSMUSG0….
Dnah14|Dnah14 41.60299 21.334681 3.4180513 6.241767 0.0000000 0.0000003 Dnah14 dynein, axonemal, heavy chain 14 240960 ENSMUSG0….
Atf3|Atf3 433.56152 -2.808330 0.7372279 -3.809311 0.0001394 0.0145805 Atf3 activating transcription factor 3 11910 ENSMUSG0….
B230208H11Rik|B230208H11Rik 48.28040 21.761939 2.7106546 8.028297 0.0000000 0.0000000 B230208H11Rik RIKEN cDNA B230208H11 gene 320273 ENSMUSG0….
LOC108167797|LOC108167797 82.55280 9.456510 2.0988898 4.505482 0.0000066 0.0009335 LOC108167797 NA NA
Fam26f|Fam26f 45.72461 21.035286 2.7780345 7.572003 0.0000000 0.0000000 Fam26f NA NA
Mical1|Mical1 109.45019 10.585487 2.6621850 3.976240 0.0000700 0.0077657 Mical1 microtubule associated monooxygenase, calponin and LIM domain containing 1 171580 ENSMUSG0….

3.2.2 Maintenance vs FoodTrained

Significant DEGs Maintenance vs FoodTrained
baseMean log2FoldChange lfcSE stat pvalue padj symbol description entrez ensembl
Gm39585|Gm39585 5797.38480 1.179803 0.2580070 4.572755 0.0000048 0.0001973 Gm39585 predicted gene, 39585 105243853 NA
Snhg6|Snhg6 522.81858 1.260214 0.4234933 2.975758 0.0029227 0.0446701 Snhg6 small nucleolar RNA host gene 6 73824 ENSMUSG0….
Trpa1|Trpa1 47.61966 -22.462571 3.6001396 -6.239361 0.0000000 0.0000001 Trpa1 transient receptor potential cation channel, subfamily A, member 1 277328 ENSMUSG0….
Rpl7|Rpl7 2486.96588 1.005240 0.3038678 3.308150 0.0009391 0.0184408 Rpl7 ribosomal protein L7 19989 ENSMUSG0….
Gm28836|Gm28836 24.57062 20.693225 3.9082681 5.294730 0.0000001 0.0000064 Gm28836 predicted gene 28836 102635593 NA
Gm32790|Gm32790 44.07127 9.563380 2.2173604 4.312957 0.0000161 0.0005913 Gm32790 predicted gene, 32790 102635458 NA
Fer1l5|Fer1l5 230.33092 -2.687277 0.8236933 -3.262473 0.0011044 0.0209149 Fer1l5 fer-1-like 5 (C. elegans) 100534273 ENSMUSG0….
Sema4c|Sema4c 1298.41081 -1.822125 0.6102459 -2.985887 0.0028276 0.0434624 Sema4c sema domain, immunoglobulin domain (Ig), transmembrane domain (TM) and short cytoplasmic domain, (semaphorin) 4C 20353 ENSMUSG0….
Vwa3b|Vwa3b 68.09876 8.717777 2.3896076 3.648204 0.0002641 0.0065801 Vwa3b von Willebrand factor A domain containing 3B 70853 ENSMUSG0….
4930439A04Rik|4930439A04Rik 239.71050 7.561383 1.1113671 6.803677 0.0000000 0.0000000 4930439A04Rik RIKEN cDNA 4930439A04 gene 78119 NA
Il1r2|Il1r2 73.49166 9.935852 2.8549805 3.480182 0.0005011 0.0112029 Il1r2 interleukin 1 receptor, type II 16178 ENSMUSG0….
Il1rl2|Il1rl2 153.74171 -24.241136 3.7848523 -6.404777 0.0000000 0.0000000 Il1rl2 interleukin 1 receptor-like 2 107527 ENSMUSG0….
LOC108167622|LOC108167622 55.54376 21.760261 2.9333442 7.418244 0.0000000 0.0000000 LOC108167622 NA NA
Gm10561|Gm10561 31.31746 -8.823060 2.8132289 -3.136275 0.0017111 0.0295377 Gm10561 predicted gene 10561 628004 NA
Clk1|Clk1 11902.12420 1.056339 0.3336618 3.165897 0.0015461 0.0272603 Clk1 CDC-like kinase 1 12747 ENSMUSG0….
Stradb|Stradb 1233.07664 1.338207 0.3402418 3.933107 0.0000839 0.0024506 Stradb STE20-related kinase adaptor beta 227154 ENSMUSG0….
Fam117b|Fam117b 1077.36822 1.740967 0.5934892 2.933444 0.0033522 0.0491723 Fam117b family with sequence similarity 117, member B 72750 ENSMUSG0….
D230017M19Rik|D230017M19Rik 353.30341 1.817985 0.5121487 3.549721 0.0003856 0.0090159 D230017M19Rik RIKEN cDNA D230017M19 gene 320933 NA
Arpc2|Arpc2 3532.07869 -1.115469 0.2391510 -4.664286 0.0000031 0.0001315 Arpc2 actin related protein 2/3 complex, subunit 2 76709 ENSMUSG0….
Cnppd1|Cnppd1 1597.07545 -1.205071 0.3232300 -3.728214 0.0001928 0.0049685 Cnppd1 cyclin Pas1/PHO80 domain containing 1 69171 ENSMUSG0….

3.2.3 Craving vs Maintenance

Significant Craving vs Maintenance
baseMean log2FoldChange lfcSE stat pvalue padj symbol description entrez ensembl
Gm39585|Gm39585 5797.38480 -1.0936780 0.2412810 -4.532799 0.0000058 0.0002749 Gm39585 predicted gene, 39585 105243853 NA
Trpa1|Trpa1 47.61966 20.0711689 3.3932805 5.914975 0.0000000 0.0000004 Trpa1 transient receptor potential cation channel, subfamily A, member 1 277328 ENSMUSG0….
Rpl7|Rpl7 2486.96588 -1.0763933 0.2841797 -3.787721 0.0001520 0.0049557 Rpl7 ribosomal protein L7 19989 ENSMUSG0….
1700001G17Rik|1700001G17Rik 704.98525 1.5933623 0.4255143 3.744557 0.0001807 0.0056755 1700001G17Rik RIKEN cDNA 1700001G17 gene 67503 ENSMUSG0….
Fer1l5|Fer1l5 230.33092 2.7365069 0.7714625 3.547168 0.0003894 0.0105650 Fer1l5 fer-1-like 5 (C. elegans) 100534273 ENSMUSG0….
Sema4c|Sema4c 1298.41081 1.9826591 0.5709143 3.472779 0.0005151 0.0131477 Sema4c sema domain, immunoglobulin domain (Ig), transmembrane domain (TM) and short cytoplasmic domain, (semaphorin) 4C 20353 ENSMUSG0….
Il1rl2|Il1rl2 153.74171 22.2736364 3.5649386 6.247972 0.0000000 0.0000001 Il1rl2 interleukin 1 receptor-like 2 107527 ENSMUSG0….
Fam126b|Fam126b 8348.97120 -0.8264248 0.2373035 -3.482565 0.0004966 0.0127789 Fam126b NA NA
Pard3bos2|Pard3bos2 48.12983 -23.7390363 3.6316177 -6.536766 0.0000000 0.0000000 Pard3bos2 par-3 family cell polarity regulator beta, opposite strand 2 105243940 ENSMUSG0….
Zdbf2|Zdbf2 3291.47772 -0.9550029 0.3014901 -3.167610 0.0015370 0.0313656 Zdbf2 zinc finger, DBF-type containing 2 73884 ENSMUSG0….
Klf7|Klf7 3114.99312 -0.6650793 0.2115799 -3.143396 0.0016700 0.0335025 Klf7 Kruppel-like factor 7 (ubiquitous) 93691 ENSMUSG0….
Map2|Map2 15113.61849 -0.6774675 0.1781094 -3.803658 0.0001426 0.0047204 Map2 microtubule-associated protein 2 17756 ENSMUSG0….
Gm31047|Gm31047 43.06422 19.1697461 3.6817645 5.206674 0.0000002 0.0000130 Gm31047 predicted gene, 31047 102633145 NA
Arpc2|Arpc2 3532.07869 0.7956339 0.2237915 3.555246 0.0003776 0.0103047 Arpc2 actin related protein 2/3 complex, subunit 2 76709 ENSMUSG0….
Cnppd1|Cnppd1 1597.07545 0.9103462 0.3024970 3.009439 0.0026173 0.0473299 Cnppd1 cyclin Pas1/PHO80 domain containing 1 69171 ENSMUSG0….
Dnajb2|Dnajb2 1487.62355 1.1910673 0.3095799 3.847366 0.0001194 0.0040449 Dnajb2 DnaJ heat shock protein family (Hsp40) member B2 56812 ENSMUSG0….
Atg9a|Atg9a 8813.08524 2.0179237 0.4649911 4.339704 0.0000143 0.0006226 Atg9a autophagy related 9A 245860 ENSMUSG0….
Resp18|Resp18 56768.10866 -0.8522294 0.2302091 -3.701980 0.0002139 0.0065238 Resp18 regulated endocrine-specific protein 18 19711 ENSMUSG0….
Chpf|Chpf 11304.97381 2.0829961 0.3845635 5.416520 0.0000001 0.0000049 Chpf chondroitin polymerizing factor 74241 ENSMUSG0….
Tmem198|Tmem198 3311.19558 1.5407766 0.3702558 4.161384 0.0000316 0.0012693 Tmem198 transmembrane protein 198 319998 ENSMUSG0….

3.3 Visualize Shrinkage Estimations

Following subsetting, I want to see how the data looks based on effect size so I have run a shrinkage model using the “apeglm” method which normalizes counts based on effect size without major changes to data structure. I have included non-shrunk models as well for comparison.

3.3.1 Non-Shrunken L2FC

3.3.1.1 Craving vs FoodTrained

Version Author Date
dc04b15 avm27 2024-07-17

3.3.1.2 Maintenance vs FoodTrained

Version Author Date
dc04b15 avm27 2024-07-17

3.3.1.3 Craving vs Maintenance

Version Author Date
dc04b15 avm27 2024-07-17

3.3.2 Shrunken L2FC

3.3.2.1 Craving vs FoodTrained

Version Author Date
dc04b15 avm27 2024-07-17

3.3.2.2 Maintenance vs FoodTrained

Version Author Date
dc04b15 avm27 2024-07-17

3.3.2.3 Craving vs Maintenance

Version Author Date
dc04b15 avm27 2024-07-17

4 Looking at Differentially Expressed Genes

We can generate a matrix to identify where our differentially expressed genes lie between each group based on an alpha value of 0.05.

## Matrix Visualizaattion
x <- vsDEGMatrix(
  data = dds, padj = 0.05, d.factor = "condition",
  type = "deseq", title = TRUE, legend = TRUE, grid = TRUE
)

Version Author Date
dc04b15 avm27 2024-07-17

5 Data Visualization

5.1 Normalized Counts and Plots of Genes of Interest

level_order <- c("FoodTrained", "Maintenance", "Craving") # this vector might be useful for other plots/analyses

normalized_counts <- counts(dds, normalized = TRUE)

5.1.1 Grin1

#### Grin1

d2 <- plotCounts(dds,
  gene = "Grin1|Grin1", intgroup = "condition",
  returnData = TRUE
)

a <- ggplot(d2, aes(x = factor(condition, level = level_order), y = count, color = condition)) +
  geom_point(position = position_jitter(w = 0.1, h = 0), size = 5) +
  scale_color_manual(values=c("#0000FF", "#FF8000", "#AD07E3")) +
  theme_bw() +
  xlab("Condition")+
  ylab("Normalized Counts")+
  ggtitle(
    label = "Grin1 Expression",
    subtitle = "Normalized Gene Expression"
  )

a

Version Author Date
21b70e3 avm27 2024-07-17
a1c090c avm27 2024-07-17
dc04b15 avm27 2024-07-17

5.1.2 Gabbr1

#### Gabbr1


d2 <- plotCounts(dds,
  gene = "Gabbr1|Gabbr1", intgroup = "condition",
  returnData = TRUE
)

a <- ggplot(d2, aes(x = factor(condition, level = level_order), y = count, color = condition)) +
  geom_point(position = position_jitter(w = 0.1, h = 0), size = 5) +
  scale_color_manual(values=c("#0000FF", "#FF8000", "#AD07E3")) +
  theme_bw() +
  xlab("Condition")+
  ylab("Normalized Counts")+
  ggtitle(
    label = "Gabbr1 Expression",
    subtitle = "Normalized Gene Expression"
  )

a

Version Author Date
21b70e3 avm27 2024-07-17
a1c090c avm27 2024-07-17
dc04b15 avm27 2024-07-17

5.1.3 Th

#### Th


d2 <- plotCounts(dds,
  gene = "Th|Th", intgroup = "condition",
  returnData = TRUE
)

a <- ggplot(d2, aes(x = factor(condition, level = level_order), y = count, color = condition)) +
  geom_point(position = position_jitter(w = 0.1, h = 0), size = 5) +
  scale_color_manual(values=c("#0000FF", "#FF8000", "#AD07E3")) +
  theme_bw() +
  xlab("Condition")+
  ylab("Normalized Counts")+
  ggtitle(
    label = "Th Expression",
    subtitle = "Normalized Gene Expression"
  )

a

Version Author Date
21b70e3 avm27 2024-07-17
a1c090c avm27 2024-07-17
dc04b15 avm27 2024-07-17

5.1.4 Ehmt2

#### Ehmt2

d2 <- plotCounts(dds,
  gene = "Ehmt2|Ehmt2", intgroup = "condition",
  returnData = TRUE
)

a <- ggplot(d2, aes(x = factor(condition, level = level_order), y = count, color = condition)) +
  geom_point(position = position_jitter(w = 0.1, h = 0), size = 5) +
  scale_color_manual(values=c("#0000FF", "#FF8000", "#AD07E3")) +
  theme_bw() +
  xlab("Condition")+
  ylab("Normalized Counts")+
  ggtitle(
    label = "Ehmt2 Expression",
    subtitle = "Normalized Gene Expression"
  )

a

Version Author Date
21b70e3 avm27 2024-07-17
a1c090c avm27 2024-07-17
dc04b15 avm27 2024-07-17

5.1.5 Setdb1

#### Setdb1

d2 <- plotCounts(dds,
  gene = "Setdb1|Setdb1", intgroup = "condition",
  returnData = TRUE
)

a <- ggplot(d2, aes(x = factor(condition, level = level_order), y = count, color = condition)) +
  geom_point(position = position_jitter(w = 0.1, h = 0), size = 5) +
  scale_color_manual(values=c("#0000FF", "#FF8000", "#AD07E3")) +
  theme_bw() +
  xlab("Condition")+
  ylab("Normalized Counts")+
  ggtitle(
    label = "Setdb1 Expression",
    subtitle = "Normalized Gene Expression"
  )

a

Version Author Date
21b70e3 avm27 2024-07-17
a1c090c avm27 2024-07-17
dc04b15 avm27 2024-07-17

5.1.6 Atf7ip

#### Atf7ip

d2 <- plotCounts(dds,
  gene = "Atf7ip|Atf7ip", intgroup = "condition",
  returnData = TRUE
)

a <- ggplot(d2, aes(x = factor(condition, level = level_order), y = count, color = condition)) +
  geom_point(position = position_jitter(w = 0.1, h = 0), size = 5) +
  scale_color_manual(values=c("#0000FF", "#FF8000", "#AD07E3")) +
  theme_bw() +
  xlab("Condition")+
  ylab("Normalized Counts")+
  ggtitle(
    label = "Atf7ip Expression",
    subtitle = "Normalized Gene Expression"
  )

a

Version Author Date
21b70e3 avm27 2024-07-17
a1c090c avm27 2024-07-17
dc04b15 avm27 2024-07-17

5.2 Principal Components Analysis (PCA)

For visualization I am converting all analyses to logarithmic. Following this I have completed a PCA for each analyses based on all genes together.

Version Author Date
21b70e3 avm27 2024-07-17
a1c090c avm27 2024-07-17
dc04b15 avm27 2024-07-17

5.3 Venn Diagrams

5.3.1 Up-Regulated Genes

myCol <- brewer.pal(3, "Pastel2")

results_sig_Table_up <- subset(results_sig, log2FoldChange > 0)
results_sig2_Table_up <- subset(results_sig2, log2FoldChange > 0)
results_sig3_Table_up <- subset(results_sig3, log2FoldChange > 0)

new_ressig1_up  <- rownames(results_sig_Table_up)
new_ressig2_up <- rownames(results_sig2_Table_up)
new_ressig3_up  <- rownames(results_sig3_Table_up)

# Generate Venn Diagram 
venn.plot.up <- venn.diagram(
  x = list(
    "Crv_vs_FT" = new_ressig1_up,
    "Maint_vs_FT" = new_ressig2_up,
    "Crv_vs_Maint" = new_ressig3_up),
    fill = myCol,
    cex = .8,
    fontface = "bold",
    fontfamily = "sans",
    cat.cex = 0.8,
    cat.fontface = "bold",
    cat.default.pos = "outer",
    cat.pos = c(0, 0, 0),
    cat.fontfamily = "sans",
    rotation = 1,
  category.names = c("Craving_vs_FoodTrained", "Maintenance_vs_FoodTrained", "Craving_vs_Maintenance"),
  filename = NULL,
  output = TRUE
)

# Save Venn Diagram
#pdf("plots/CIVSA_venn_diagram_Up.pdf")
grid.draw(venn.plot.up)

Version Author Date
dc04b15 avm27 2024-07-17
#dev.off()

# Display Venn Diagram in RStudio
#grid.draw(venn.plot.up)

5.3.2 Down-Regulated Genes

results_sig_Table_down <- subset(results_sig, log2FoldChange < 0)
results_sig2_Table_down <- subset(results_sig2, log2FoldChange < 0)
results_sig3_Table_down <- subset(results_sig3, log2FoldChange < 0)

new_ressig1_down  <- rownames(results_sig_Table_down)
new_ressig2_down  <- rownames(results_sig2_Table_down)
new_ressig3_down  <- rownames(results_sig3_Table_down)


# Generate Venn Diagram UP
venn.plot.down <- venn.diagram(
  x = list(
    "Crv_vs_FT" = new_ressig1_down,
    "Maint_vs_FT" = new_ressig2_down,
    "Crv_vs_Maint" = new_ressig3_down),
    fill = myCol,
    cex = .8,
    fontface = "bold",
    fontfamily = "sans",
    cat.cex = 0.8,
    cat.fontface = "bold",
    cat.default.pos = "outer",
    cat.pos = c(0, 0, 0),
    cat.fontfamily = "sans",
    rotation = 1,
  category.names = c("Craving_vs_FoodTrained", "Maintenance_vs_FoodTrained", "Craving_vs_Maintenance"),
  filename = NULL,
  output = TRUE
)

# Save Venn Diagram
#pdf("plots/CIVSA_venn_diagram_Down.pdf")
#grid.draw(venn.plot.down)
#dev.off()

# Display Venn Diagram in RStudio
grid.draw(venn.plot.down)

Version Author Date
dc04b15 avm27 2024-07-17

5.4 Heatmaps of Differentially Expressed Genes

ressig1 <- as.data.frame(results_sig)
ressig2 <- as.data.frame(results_sig2)
ressig3 <- as.data.frame(results_sig3)

ressig1$Comparison <- "Craving vs FoodTrained"
ressig2$Comparison <- "Maintenance vs FoodTrained"
ressig3$Comparison <- "Craving vs Maintenance"

ressig1_highlfc_pos <- ressig1[ressig1$log2FoldChange >= 1.5,]
ressig1_highlfc_neg <- ressig1[ressig1$log2FoldChange <= -1.5,]
ressig2_highlfc_pos <- ressig2[ressig2$log2FoldChange >= 1.5,]
ressig2_highlfc_neg <- ressig2[ressig2$log2FoldChange <= -1.5,]
ressig3_highlfc_pos <- ressig3[ressig3$log2FoldChange >= 1.5,]
ressig3_highlfc_neg <- ressig3[ressig3$log2FoldChange <= -1.5,]

temp_sigdegshm <- rbind(ressig1_highlfc_pos, ressig1_highlfc_neg)
temp_sigdegshm2 <- rbind(temp_sigdegshm, ressig2_highlfc_pos)
temp_sigdegshm3 <- rbind(temp_sigdegshm2, ressig2_highlfc_neg)
temp_sigdegshm4 <- rbind(temp_sigdegshm3, ressig3_highlfc_pos)
temp_sigdegshm5 <- rbind(temp_sigdegshm4, ressig3_highlfc_neg)

temp_sigdegshm6 <- temp_sigdegshm5[temp_sigdegshm5$lfcSE <=1,]

dup_rows <- duplicated(temp_sigdegshm6$symbol)
data_dedup <- temp_sigdegshm6[!dup_rows, ]

rows <- rownames(data_dedup)

select <- order(rowMeans(counts(dds,normalized=TRUE)),
                decreasing=TRUE)#[1:1500]

mat <- assay(ddsMat_rlog)[select,]

data_mod <- mat[rownames(mat) %in% rows, ] 

matrownames_data_mod <- data.frame(rownames(data_mod))

matrownames_data_mod <- matrownames_data_mod %>% separate(rownames.data_mod., c("id", "symbol"),sep ="([|])")

samples <- as.tibble(read.csv("data/sample_info_CIVSA.csv", sep = ","))

annot_col <- samples %>%
  column_to_rownames('samplename') %>%
  as.data.frame()

ann_colors = list(
  condition = as.factor(c(FoodTrained = "#0000FF", Maintenance = "#AD07E3", Craving = "#FF8000")),
  Group = as.factor(c(FoodTrained = "#0000FF", Maintenance = "#AD07E3", Craving = "#FF8000")))

5.4.1 Unsupervised Clustering

pheatmap::pheatmap(
  mat = data_mod,
  color = colorRampPalette(rev(brewer.pal(11, "RdYlBu")))(155),
  scale = "row",
  annotation_colors = ann_colors,
  annotation_col = annot_col,
  fontsize = 10, 
  show_colnames = T,
  cluster_rows = T,
  cluster_cols = T,
  annotation_legend = T,
  legend = T,
  show_rownames = F,
  clustering_distance_rows = "correlation",
  clustering_method = "ward.D",
  main = "Heatmap of Significant DEGs -- Unsupervised Clustering",
  cluster_row_slices = TRUE,
  show_row_dend = FALSE,
  cutree_rows = 4,
  cutree_cols = 3,
  treeheight_row = 10,
  treeheight_col = 10,
  #legend_labels = c("FoodTrained", "Maintenance", "Craving")
)

Version Author Date
21b70e3 avm27 2024-07-17
a1c090c avm27 2024-07-17
dc04b15 avm27 2024-07-17

5.4.2 Supervised Clustering

pheatmap::pheatmap(
  mat = data_mod,
  color = colorRampPalette(rev(brewer.pal(11, "RdYlBu")))(155),
  scale = "row",
  annotation_col = annot_col,
  annotation_colors = ann_colors,
  fontsize = 10,
  show_colnames = T,
  cluster_rows = T,
  cluster_cols = F,
  annotation_legend = T,
  legend = T,
  show_rownames = F,
  clustering_distance_rows = "correlation",
  clustering_method = "ward.D",
  main = "Heatmap of Significant DEGs -- Supervised Clustering",
  cluster_row_slices = TRUE,
  show_row_dend = FALSE,
  cutree_rows = 4,
  cutree_cols = 3,
  gaps_col = c(3, 6),
  treeheight_row = 10,
  treeheight_col = 10,
  #legend_labels = c("FoodTrained", "Maintenance", "Craving"),
)

Version Author Date
21b70e3 avm27 2024-07-17
a1c090c avm27 2024-07-17
dc04b15 avm27 2024-07-17

5.5 Volcano Plots

5.5.1 Visualize Results: Volcano Plots (Non-Interactive)

Volcano plots are generated by gathering the FDR corrected p-values from each analysis. The adjusted p-values undergo a -log10 transformation to generate these FDR corrected values. Any rows containing “NA” are ommitted from analysis. Data points are colored based on increasing or decreasing value and plotted using ggplot2. Adjusted p-value cutoff is 1.3 following log transformation.(-log10 0.05 = 1.3)

Version Author Date
a1c090c avm27 2024-07-17

5.5.2 Visualize Results: Volcano Plots (Interactive)

Interactive Volcano plots are generated by Using the same method as above, but are plotted so that genes may be identified and log2FC can be seen from each gene. Plots can be used to customize which genes appear based on what you would like to see.


sessionInfo()
R version 4.2.1 (2022-06-23)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 22.04.4 LTS

Matrix products: default
BLAS:   /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.10.0
LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.10.0

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] grid      stats4    stats     graphics  grDevices utils     datasets 
[8] methods   base     

other attached packages:
 [1] plotly_4.10.2                            
 [2] gridExtra_2.3                            
 [3] ggrepel_0.9.3                            
 [4] BiocParallel_1.32.6                      
 [5] conflicted_1.2.0                         
 [6] VennDiagram_1.7.3                        
 [7] futile.logger_1.4.3                      
 [8] regionReport_1.32.0                      
 [9] kableExtra_1.3.4                         
[10] vidger_1.18.0                            
[11] tables_0.9.17                            
[12] knitr_1.43                               
[13] Glimma_2.8.0                             
[14] Mus.musculus_1.3.1                       
[15] TxDb.Mmusculus.UCSC.mm10.knownGene_3.10.0
[16] OrganismDbi_1.40.0                       
[17] GenomicFeatures_1.50.4                   
[18] gplots_3.1.3                             
[19] lubridate_1.9.2                          
[20] forcats_1.0.0                            
[21] stringr_1.5.1                            
[22] purrr_1.0.2                              
[23] readr_2.1.4                              
[24] tibble_3.2.1                             
[25] tidyverse_2.0.0                          
[26] statmod_1.5.0                            
[27] tweeDEseqCountData_1.36.0                
[28] edgeR_3.40.2                             
[29] limma_3.54.2                             
[30] GOSemSim_2.24.0                          
[31] apeglm_1.20.0                            
[32] ggpubr_0.6.0                             
[33] SPIA_2.50.0                              
[34] KEGGgraph_1.58.3                         
[35] ggnewscale_0.4.9                         
[36] enrichplot_1.18.4                        
[37] tidyr_1.3.1                              
[38] WGCNA_1.72-1                             
[39] fastcluster_1.2.3                        
[40] dynamicTreeCut_1.63-1                    
[41] pathview_1.38.0                          
[42] gage_2.48.0                              
[43] dplyr_1.1.4                              
[44] topGO_2.50.0                             
[45] SparseM_1.81                             
[46] graph_1.76.0                             
[47] GO.db_3.16.0                             
[48] RColorBrewer_1.1-3                       
[49] genefilter_1.80.3                        
[50] pheatmap_1.0.12                          
[51] ggsci_3.0.0                              
[52] tximport_1.26.1                          
[53] org.Mm.eg.db_3.16.0                      
[54] AnnotationDbi_1.60.2                     
[55] DOSE_3.24.2                              
[56] ReactomePA_1.42.0                        
[57] biomaRt_2.54.1                           
[58] clusterProfiler_4.6.2                    
[59] ggplot2_3.5.0.9000                       
[60] DESeq2_1.38.3                            
[61] SummarizedExperiment_1.28.0              
[62] Biobase_2.58.0                           
[63] MatrixGenerics_1.10.0                    
[64] matrixStats_1.0.0                        
[65] GenomicRanges_1.50.2                     
[66] GenomeInfoDb_1.34.9                      
[67] IRanges_2.32.0                           
[68] S4Vectors_0.36.2                         
[69] BiocGenerics_0.44.0                      
[70] gprofiler2_0.2.2                         
[71] workflowr_1.7.0                          

loaded via a namespace (and not attached):
  [1] Hmisc_5.1-0              svglite_2.1.1            ps_1.7.5                
  [4] Rsamtools_2.14.0         foreach_1.5.2            rprojroot_2.0.3         
  [7] crayon_1.5.2             MASS_7.3-60              nlme_3.1-162            
 [10] backports_1.4.1          impute_1.72.3            rlang_1.1.3             
 [13] XVector_0.38.0           HDO.db_0.99.1            callr_3.7.3             
 [16] filelock_1.0.2           rjson_0.2.21             bit64_4.0.5             
 [19] glue_1.7.0               rngtools_1.5.2           parallel_4.2.1          
 [22] processx_3.8.2           DEFormats_1.26.0         tidyselect_1.2.1        
 [25] XML_3.99-0.14            GenomicAlignments_1.34.1 xtable_1.8-4            
 [28] magrittr_2.0.3           evaluate_0.21            bibtex_0.5.1            
 [31] cli_3.6.2                zlibbioc_1.44.0          doRNG_1.8.6             
 [34] rstudioapi_0.14          whisker_0.4.1            bslib_0.5.0             
 [37] rpart_4.1.19             derfinderHelper_1.32.0   fastmatch_1.1-3         
 [40] BiocStyle_2.26.0         lambda.r_1.2.4           treeio_1.22.0           
 [43] xfun_0.39                gson_0.1.0               cluster_2.1.4           
 [46] caTools_1.18.2           tidygraph_1.2.3          KEGGREST_1.38.0         
 [49] ape_5.7-1                Biostrings_2.66.0        png_0.1-8               
 [52] reshape_0.8.9            withr_3.0.0              bitops_1.0-7            
 [55] ggforce_0.4.1            RBGL_1.74.0              plyr_1.8.8              
 [58] coda_0.19-4              bumphunter_1.40.0        pillar_1.9.0            
 [61] cachem_1.0.8             fs_1.6.2                 graphite_1.44.0         
 [64] ellipsis_0.3.2           vctrs_0.6.5              generics_0.1.3          
 [67] tools_4.2.1              foreign_0.8-84           munsell_0.5.1           
 [70] tweenr_2.0.2             fgsea_1.24.0             DelayedArray_0.24.0     
 [73] fastmap_1.1.1            compiler_4.2.1           abind_1.4-5             
 [76] httpuv_1.6.11            rtracklayer_1.58.0       GenomeInfoDbData_1.2.9  
 [79] lattice_0.21-8           utf8_1.2.4               later_1.3.1             
 [82] BiocFileCache_2.6.1      jsonlite_1.8.7           GGally_2.1.2            
 [85] scales_1.3.0             tidytree_0.4.2           carData_3.0-5           
 [88] lazyeval_0.2.2           promises_1.2.0.1         car_3.1-2               
 [91] doParallel_1.0.17        checkmate_2.2.0          rmarkdown_2.23          
 [94] cowplot_1.1.1            webshot_0.5.5            downloader_0.4          
 [97] BSgenome_1.66.3          igraph_1.5.0             survival_3.5-5          
[100] numDeriv_2016.8-1.1      yaml_2.3.7               systemfonts_1.0.4       
[103] htmltools_0.5.5          memoise_2.0.1            VariantAnnotation_1.44.1
[106] BiocIO_1.8.0             locfit_1.5-9.8           graphlayouts_1.0.0      
[109] viridisLite_0.4.2        digest_0.6.32            rappdirs_0.3.3          
[112] knitrBootstrap_1.0.2     futile.options_1.0.1     emdbook_1.3.13          
[115] RSQLite_2.3.1            yulab.utils_0.0.6        derfinder_1.32.0        
[118] data.table_1.14.8        blob_1.2.4               preprocessCore_1.60.2   
[121] labeling_0.4.3           splines_4.2.1            Formula_1.2-5           
[124] RCurl_1.98-1.12          broom_1.0.5              hms_1.1.3               
[127] colorspace_2.1-0         base64enc_0.1-3          BiocManager_1.30.21     
[130] aplot_0.1.10             nnet_7.3-19              sass_0.4.6              
[133] Rcpp_1.0.10              mvtnorm_1.2-2            fansi_1.0.6             
[136] tzdb_0.4.0               R6_2.5.1                 lifecycle_1.0.4         
[139] formatR_1.14             curl_5.0.1               ggsignif_0.6.4          
[142] jquerylib_0.1.4          Matrix_1.5-4.1           qvalue_2.30.0           
[145] org.Hs.eg.db_3.16.0      iterators_1.0.14         RefManageR_1.4.0        
[148] htmlwidgets_1.6.2        markdown_1.7             polyclip_1.10-4         
[151] crosstalk_1.2.0          shadowtext_0.1.2         timechange_0.2.0        
[154] gridGraphics_0.5-1       reactome.db_1.82.0       rvest_1.0.3             
[157] htmlTable_2.4.1          patchwork_1.2.0.9000     bdsmatrix_1.3-6         
[160] codetools_0.2-19         gtools_3.9.4             getPass_0.2-2           
[163] prettyunits_1.1.1        dbplyr_2.3.2             gtable_0.3.5            
[166] DBI_1.1.3                git2r_0.32.0             ggfun_0.1.1             
[169] httr_1.4.6               highr_0.10               KernSmooth_2.23-21      
[172] stringi_1.8.4            progress_1.2.2           reshape2_1.4.4          
[175] farver_2.1.2             annotate_1.76.0          viridis_0.6.3           
[178] Rgraphviz_2.42.0         ggtree_3.6.2             xml2_1.3.4              
[181] bbmle_1.0.25             restfulr_0.0.15          geneplotter_1.76.0      
[184] ggplotify_0.1.1          bit_4.0.5                scatterpie_0.2.1        
[187] ggraph_2.1.0             pkgconfig_2.0.3          rstatix_0.7.2           
[190] GenomicFiles_1.34.0